Fractional Bloch oscillations in photonic lattices 
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Bloch oscillations (BO), i.e., the oscillatory motion of a quantum particle in a periodic potential 
driven by a constant force, constitute one of the most striking and oldest predictions of coherent 
quantum transport in periodic lattices. In natural crystals, BO have never been observed because 
of dephasing effects. Solely with the advent of semiconductor super lattices and ultracold atoms BO 
have been observed for matter waves [iJB]. In their essence, BO are a wave phenomenon. As such, 
they are found for optical [SHS] and acoustic [10] waves as well. When interactions between particles 
compete with their mobility, novel dynamical behaviour can arise where particles form bound states 
pT| and co-tunnel through the lattice [12] . While particle interaction has been generally associated 
to BO damping 13H15j. for few strongly-interacting particles it was predicted that bound states 
undergo fractional BO at a frequency twice (or multiple) that of single-particle BO [16-18 . The 
observation of fractional BO is challenging in condensed-matter systems and up to now has not been 
achieved even in model systems. Here we report on the first observation of fractional BO using a 
photonic lattice as a model system of a few-particle extended Bose-Hubbard Hamiltonian. These 
results pave the way to the visualization of other intriguing phenomena involving few correlated 
particles, such as the interplay between particle interaction and Anderson localization 1191 120) . 
quantum billiards and transition to quantum chaos [21], anyonic BO [22 , and correlated barrier 
tunneling and particle dissociation [23] . 



Ultracold quantum gases in optical lattices have pro- 
vided in the past decade a powerful laboratory tool to 
simulate Hubbard models of condensed-matter physics 
[24] . For few interacting particles, experiments with ul- 
tracold atoms have so far demonstrated the existence of 
bound particle states and correlated tunneling phenom- 
ena |lTJ[T2], however fractional BO of atomic pairs have 
not been observed yet. A photonic simulator of the Hub- 
bard model can map the dynamics of two correlated par- 
ticles hopping on a one dimensional lattice into the mo- 
tion of a single particle in a two-dimensional lattice with 
engineered defects [2TJ [25], thereby offering the possi- 
bility to easily visualize the effect under simulation, but 
requiring a high level of control of the photonic structure. 

To study correlated BO, we consider an extended Bose- 
Hubbard (EBH) model [^G', 'ST describing strongly inter- 
acting bosons in the lowest Bloch band of a one dimen- 
sional lattice driven by an external force. This model 
is more accurate than the standard Bose-Hubbard (BH) 
model as it accounts for higher-order processes whose 
magnitude is comparable with the one of second-order 
tunneling [27] . The Hamiltonian of the system is given 
hy H = Hebh + Hp^ where Hp = Fd^^ Ifii describes 
the effect of the external constant force F {d is the lattice 
period) and Hebh = Hbh -\- Hs -\- H4 -\- is the EBH 
Hamiltonian. Hebh includes the standard BH Hamilto- 
nian (with h = 1) 

Hbh = y (^^ - 1) - ^ {^l^i+i + H.c.) (1) 
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and three additional terms [26] , whose explicit expression 
is given in the Supplementary information. In Eq.(l), 
Uq is the on-site energy interaction strength {Uq > 
for a repulsive interaction); J is related to the single- 
particle hopping rate between adjacent lattice sites; e <C 
1 is the lattice attenuation factor (see [26 ); dj and 
di are the bosonic creation and annihilation operators 
and hi = d\di is the particle number operator at lat- 
tice site /. The standard BH term Hbh describes on- 
site particle-particle interaction and single-particle tun- 
neling between adjacent lattice sites [see Figs. 1(a) and 
(b)]. The additional terms H^,^ H4 and H^ account for 
nearest-neighbor hopping of single atoms conditioned by 
the on-site occupation number, nearest-neighbor inter- 
actions, and nearest-neighbor hopping of bosonic pairs, 
respectively [see Figs. 1(c) and (d)]. 

In the absence of particle interaction, each particle 
undergoes independent BO on the lattice. The single- 
particle amplitude probabilities An{t) to find the particle 
at lattice site n evolve according to the coupled equations 

dA 

i-^ = -^(An+i + An-i) + FdnAn (2) 

where n = eJ/2 is the single-particle hopping rate. The 
single-particle energy spectrum is given by an equally- 
spaced Wannier- Stark ladder, which results in a periodic 
dynamics at the BO frequency ojb = Fd. For single-site 
excitation, i.e. An(0) = ^n,0 7 BO appear as a breathing 
(rather than an oscillatory) motion at the frequency uob- 
Photonic demonstrations of such BO breathing dynamics 
were previously reported in Refs.fS', '6^, 'S . 

Let us now consider the dynamics of two interact- 
ing particles. The state vector \ip{t)) of the sys- 
tem can be expanded in Fock space as |'0(t)) = 
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FIG. 1. Dynamics of two interacting bosons in the extended Bose-Hubbard lattice. Two-particle dynamics in the 
EBH model: (a) Single-particle tunneling; (b) on-site particle interaction; (c) nearest-neighbor-site particle interaction and 
conditional single-particle tunneling; (d) nearest- neighbor hopping of particle pairs, (e) Fock-space representation of the two- 
particle dynamics; the upper-right inset represents the two possible pathways for two-particle hopping to the nearest-neighbor 
site: (a) is the second-order pair tunneling, while is the direct two-particle tunneling. 



{l/V^)J2n,m'^n,m{t)alal,\0), where Cn,m{t) is the am- 
plitude probability to find one particle at site n and the 
other particle at site m, with Cn,m = Cm,n for bosons. The 
evolution equations for the amplitude probabilities (see 
the Supplementary information) show that the dynam- 
ics of the two correlated bosons on the one-dimensional 
lattice can be mapped into the motion of a single par- 
ticle hopping on the two-dimensional square lattice of 
Fig. 1(e). The square lattice shows both energy site 
defects and hopping rate corrections on the three di- 
agonals m = n,n ± 1. Far from the main diagonal, 
the square lattice is homogeneous with a uniform hop- 
ping rate = eJ/2, which is the single-particle hopping 
rate of the original problem [see Fig. 1(a)]. The site en- 
ergy defects Uq and 2e^Uo on the main (m = n) and 
nearest-neighbor (m = nil) diagonals account for on- 
site and nearest-neighbor-site particle repulsion, respec- 
tively [Figs. 1(b) and (c)]. The correction of hopping 
rates k,i = k, -\- k,' (with k,' = —Uqc^^'^) on the main 
diagonal accounts for nearest-neighbor hopping of sin- 
gle atoms conditioned by the on-site occupation num- 
ber [Fig. 1(c)]. Finally, the cross coupling p = —2Uoe'^ 
on the main diagonal n = m describes nearest-neighbor 
hopping of bosonic pairs, i.e. direct two-particle tunnel- 
ing [Fig. 1(d)]. For BO of two strongly-correlated par- 
ticles (J ~ t/o), only this last term is important, com- 
peting with second-order pair tunneling. In fact, let us 
assume as an initial condition Cn,m(0) = An{0)Sn,m^ i-e. 
that the two particles are initially placed on the same 
site, with a probability [^^(O)^ to find both particles 



at site n. Then an asymptotic analysis of the underly- 
ing equations shows that the two particles form a bound 
state, i.e. they co-tunnel along the lattice, undergo- 
ing BO at the frequency 2ujb = 2Fd twice the single- 
particle BO frequency uob- In Fock space, this means 
that the excitation remains confined along the main di- 
agonal m = n of the square lattice of Fig. 1(e), and pe- 
riodic revivals at frequency 2ujb are expected. Precisely, 
as shown in the Supplementary information, the bound- 
particle occupation amplitudes An{t) ~ Cn,n(^) at the 
various lattice sites evolve according to Eqs.(2), but with 
F replaced by 2F and n replaced by i^eff- The effec- 
tive hopping rate of the bound particle state is given 
by keff = -2i^^/Uo + p = -e^jy{2Uo) - 2Uoe^ , which 
includes second-order tunneling {—2hz'^/Uo) and direct 
two-boson tunneling (p) processes. Such two contribu- 
tions always add with the same sign, for both repulsive 
{Uo > 0) and attractive {Uq < 0) interactions, and gen- 
erally are of the same order of magnitude. The bound 
particle state thus behaves like a single particle hopping 
on the one-dimensional lattice, but with a driving force 
which is doubled as compared to that of a single particle 
[T6| [TT] and with a modified hopping rate i^eff- Note 
that the effect of direct two- atom tunneling on the BO 
dynamics, not considered in previous works [HI [181 [25] , 
is to increase the amplitude of the breathing motion (ow- 
ing to the correction of i^eff)^ whereas the Wannier-Stark 
energy spectrum is not altered, i.e. cj^ is not affected. 

In our experiments, single-particle hopping dynamics 
in the square lattice of Fig. 1(e) has been simulated by 
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FIG. 2. The photonic simulator of correlated BO. (a) Sketch of the waveguide array structure; red-coloured waveguides 
have a different refractive index change with respect to the others to implement the on-site particles-interaction defect Uq. (b) 
Section of the fabricated array, imaged with an optical microscope. 



discretized spatial light transport in an engineered two- 
dimensional square lattice of evanescently-coupled opti- 
cal waveguides [28^. This lattice has been fabricated in 
a fused silica substrate by direct waveguide writing with 
femtosecond lasers [29^ "33 , taking advantage of the three- 
dimensional capabilities of this technology. The spatial 
propagation of the light intensity in the generic (m, n)- 
th lattice site maps the temporal evolution of the two- 
particle probability distribution |c^,n(^)P in Fock space 
[28[ [32] . In the photonic simulator of HebHi the role 
of the hopping coefficients n and p is played by the 
evanescent coupling constants between first- and second- 
neighboring waveguides {R and p), respectively. The for- 
mer is in general greater than the latter and the relative 
weight can be tuned by properly engineering the lattice 
parameter d (see Supplementary Information). The on- 
site particle-interaction defect [/q is realized by fabricat- 
ing the main diagonal waveguides with a different refrac- 
tive index change (this is achieved by slightly modifying 
the writing speed of the waveguides) and thus with the 
propagation constant P' slightly detuned with respect to 
that of all the others (/3). The positive or negative sign 
of the relative detuning A/3 = P — j3' determines the 
repulsive or attractive nature of the interaction, respec- 
tively. Lastly, we did not add any special feature in the 
nearest-neighboring diagonals, because, as already dis- 
cussed above, the nearest-site effects described in Fig. 1(c) 
are negligible with respect to the other terms of the EBH 
model. The implementation of the driving force term Hp 
is obtained by a constant bending of the waveguide axes 
(radius of curvature F ex 1/i^) in the plane deter- 
mined by the lattice main diagonal direction Y and the 
light propagation direction z [8 . In Fig. 2 (a) the whole 
lattice structure is depicted. Note that the array cross- 
section is rotated in such a way that the bending plane is 
horizontal in the laboratory reference frame [45° rotation 
with respect to Fig. 1(e)]. This simplifies the fabrication 
process and the imaging of the light propagating in the 
waveguides; in this way, although being a 3D structure, 



each waveguide remains on a specific horizontal plane. 

In a first experiment we checked that our photonic 
simulator can reliably mimic the dynamics governed by 
Hebh alone, namely without the action of the exter- 
nal constant force. For this purpose we fabricated 225 
straight waveguides {R = oo) arranged in a 15x15 square 
lattice, with spacing d = Idjam and a length L = 2.5cm. 
For such a lattice spacing the coupling coefficients are 
R = 0.95cm~^ and p = 0.3cm~^. To simulate particle 
interaction, we fabricated the main diagonal waveguides 
with a propagation constant detuning A/3 = — 4cm~^, 
therefore representing an attractive interaction between 
the particles. A microscope image of the array section is 
reported in Fig. 2(b). We coupled a 633nm He-Ne laser 
radiation into the central waveguide, belonging to the 
main diagonal subset and thus corresponding to an ini- 
tial condition where the two particles occupy the same 
lattice site. As previously discussed light is expected to 
remain on the detuned diagonal, representing the fact 
that the correlated particles will hop together to the next 
site and will never separate. To study the evolution of 
the system we monitored the light propagation along the 
array by imaging from above the fluorescence light at 
650 nm emitted from the waveguides and arising from 
colour centres created in the femtosecond laser writing 
process. Since the colour centres are formed exclusively 
inside the waveguides, this technique yields a high signal- 
to-noise ratio and permits a direct imaging of the light 
propagation in the waveguide arrays [32l [33] . The result 
obtained in this first experiment is visible in Fig. 3(a). 
The measured distribution of light in the photonic model 
implementing both contributions to the bound state hop- 
ping rate, i.e. two-particle tunneling and second-order 
pair tunneling, is compared with the corresponding nu- 
merical simulation [Fig. 3(b)]. The excellent agreement 
demonstrates that our photonic model can indeed take 
into account both contributions and thus correctly im- 
plements the EBH model. The role of the direct two- 
particle tunneling [Fig. 1(d)], neglected in the standard 
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FIG. 3. Photonic simulation of the delocalization dy- 
namics of two interacting particles in a crystal, (a) 

Experimental simulation in a waveguide lattice implementing 
direct two-particle tunneling and second-order pair tunneling 
with a similar weight, (b) Numerical simulation including 
both contributions, (c) Numerical simulation including only 
second-order pair tunneling. 



BH model, can be quite relevant as shown by the different 
light distributions obtained from numerical simulations 
of the same structure but including only the second-order 
tunneling [Fig. 3(c)]. 

In a second experiment, we proceeded to the observa- 
tion of fractional BO. To this purpose we fabricated an- 
other 15x15 square waveguide lattice, with the same lat- 
tice spacing and detuning on the diagonal as in the previ- 
ous structure. To mimic the external force F, the waveg- 
uides are now circularly-bent with a radius R = 400cm, 
as depicted in Fig. 2 (a). The array length is L = 8.5cm. 
As in the previous experiment, we injected the probe light 
into the central waveguide and we imaged its propagation 
along the waveguide lattice diagonal (Y), where light is 
confined. Figure 4(a) shows the oscillatory behaviour of 
the light dynamics along the waveguide array. The accu- 
racy of the fabricated photonic simulator is confirmed by 
the very good agreement with the numerical simulations 
of light propagation in the designed structure [Fig. 4(a)]. 
Initially, light spreads into several waveguides, until it 
reaches a maximum of the breathing amplitude (around 
Lmax = 3.25cm). Then light refocuses into the cen- 
tral waveguide {Lfoc = 6.5cm) and the breathing starts 
again. This periodic behavior represents the correlated 
BO of two interacting particles, watched in Fock space. 
To better analyze this result we cut the sample at Lmax 
and Lfoc and we imaged onto a CCD camera the light 
distribution exiting the output facet of the array (Fig. 5). 
Figure 5(a) clearly shows that, even in correspondence of 




FIG. 4. Photonic simulation of Bloch oscillations for 
correlated and single particles, a) Experimental and nu- 
merical simulation of the light distribution in the waveguide 
lattice diagonal, representing in Fock space the BO for two 
interacting particles . (b) Experimental and numerical simu- 
lation of the light distribution corresponding to single particle 
BO; the position where the BO amplitude is at maximum cor- 
responds to the refocusing condition for the two-particle BO 
Lfoc — 6.5cm. 



the breath maximum amplitude, light remains confined 
in the waveguide lattice detuned diagonal. This fact is a 
further demonstration that two interacting particles, ini- 
tially in the same lattice site, form a bound state and hop 
together to neighboring sites. Figure 5(b) shows that at 
the refocusing distance, all the light is indeed confined 
in the central waveguide, which was the one originally 
excited (indicated by a white arrow in both panels). 

The 'fractional' nature of BO for correlated particles, 
namely frequency doubling with respect to the case of 
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FIG. 5. Output light distribution at two propagation 
lengths. Measured light distribution in the waveguide lat- 
tice section [corresponding to Fig. 4(a)] at positions (a) Lmax, 
where the fractional BO amplitude is at maximum and (b) 
Lfoc where light refocuses to the initial waveguide (indicated 
by a white arrow in both panels). 

a single particle hopping in the same lattice under the 
action of the same force F, is shown in Fig. 4(b). We 
fabricated a planar ID array of 23 identical waveguides, 
equally spaced of d' = d = 19/im, and uniformly bent in 
the array plane with a constant radius of curvature R\ 
Such a structure reliably simulates the dynamics of sin- 
gle particle BO [8 . In order to implement the same force 
on the single particle we have to project the curvature 
on the diagonal for the two-particles structure on one of 
the main axis of the square lattice: this yields a curva- 
ture radius for the linear array 

We imaged the light propagation in this waveguide ar- 
ray for a length L' = L = 8.5cm [Fig. 4(b)]. Also in 
this case we observe a breathing propagation mode cor- 
responding to the single particle BO, but the position of 
the maximum breath amplitude L'^^^ corresponds to the 
refocusing position Lfoc for two-particles BO [Fig. 4(a)]. 
This demonstrates that the frequency of a two-particles 
BO is twice that for a single-particle. It may be worth 
noting that the breathing amplitude of a single particle 
BO is larger than that for the two-particle BO. This is 
due to a higher hopping rate in the former case, k > keff^ 
and to the double force experienced by the two strongly 
correlated particles with respect to the single one. 

In conclusion, we have experimentally observed the 
photonic quantum- analogy of fractional Bloch Oscilla- 
tions of two strongly correlated particles. We have shown 
that our photonic simulator is capable of implementing 
two-particle tunneling in addition to second-order pair 
tunneling, thus enabling the simulation of the extended 
Bose-Hubbard model. The simple visualization of the 
phenomenon under study and the high control on the 
fabricated structure, and thus of the implemented model, 
make this approach a very powerful tool to investigate 
other exotic phenomena associated to the formation of 
particle bound states which are of difficult access in the 
matter. 
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METHODS 

Fabrication of the waveguide arrays. Waveguide ar- 
rays have been fabricated by femtosecond laser waveguide 
writing in fused silica samples. The second harmonic of a 
HighQLaser femtoREGEN Yb-based amplified laser system 
has been used, consisting in 400- fs pulses at 520-nm wave- 
length with repetition rates up to 960 kHz. Actual writing 
conditions consisted in 300-nJ pulses delivered at a repetition 
rate of 20 kHz, reduced using the internal pulse-picker. Opti- 
mal range for the writing speed was identified as 8 — 14mm/s, 
yielding propagation losses of about 0.5dB/cm. In order to 
implement the detuning in the lattice diagonal we modulated 
the writing speed, therefore the waveguides on the diagonal 
were fabricated at a speed of 9mm/ s, while all the others were 
written at 14mm/s. The variation in writing speed modifies 
the propagation constants of the guided modes, producing the 
designed detuning A/3 = — 4cm~^, but has a negligible effect 
on the coupling rate between the waveguides. End-faces of 
the substrates were polished after writing to improve light 
launching and output imaging. 

Characterization of the light distribution in the ar- 
rays. Femtosecond laser writing in fused silica creates color 
centers that provide fluorescent emission at about 650 nm 
when light at 633 nm is propagated in the waveguide. Top- 
view imaging of the fluorescence signal is employed to visual- 
ize and quantitatively estimate the light distribution along the 
waveguide array, rejecting the background light by a notch fil- 
ter at 633 nm. In order to achieve a high resolution in imaging 
the light in the waveguide array all along its length, several 
images have been acquired and then stitched together. Prop- 
agation losses in the waveguides have been compensated by 
renormalizing the intensity levels in the acquired images. 
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